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ABSTRACT 

We present the first kinematic analysis of the far outer halo globular cluster (GC) population in the Local 
Group galaxy M31. Our sample contains 53 objects with projected radii of ^20-130 kpc, of which 44 have 
no previous spectroscopic information. GCs with projected radii > 30 kpc are found to exhibit net rotation 
around the minor axis of M31, in the same sense as the inner GCs, albeit with a smaller amplitude of 79+19 
km/s. The rotation-corrected velocity dispersion of the full halo GC sample is 106+12 km/s, which we observe 
to decrease with increasing projected radius. We find compelling evidence for kinematic-coherence amongst 
GCs which project on top of halo substructure, including a clear signature of infall for GCs lying along the 
North- West stream. Using the tracer mass estimator, we estimate the dynamical mass of M31 within 200 kpc 
to be Mm3i = (1 .2- 1 .5) ± 0.2 x IO'^Mq. This value is highly dependent on the chosen model and assumptions 
within. 

Subject headings: Local Group — galaxies: individual (M31) — galaxies: kinematics and dynamics — galax- 
ies: halos — globular clusters: general 



1. INTRODUCTION 

Globular cluster (GC) systems contai n important clue s 
about the assembly history of galaxies (e.g.. West et al. 2004). 
Their kinematics are especially important as dif ferent forma- 
tion channels lead to distinc t predictions (e.g. iForbes et al.l 
ll997l:lBrodie & Straderil2006l) . Moreover, GC kinematics can 
also be used to model the shape of the gravitation al potential 
and c onstrain the total mass of the host galaxy (e.g. lCote et al.l 
l200l . 

Loc ated at a distance of ^780 kpc (iMcConnachie et al] 
12005b . the Local Group galaxy M3 1 provides an excellent op- 
portunity to study a rich GC system in unparalleled detail. 
It has more than 500 confirmed members listed in the Re- 
vised Bologna Catalogue (RBC, Galleti et al. 2004), most of 
which lie within a projected radius (Rproj) of 3 kpc. In recent 
years , state-of-the-art w i de field surveys (e.g. [Ferguson et al.l 
I200I libataetallllOOTl: IMcConnachie et alj|2009h have en- 
abled searches for GCs in M3rs far outer halo. This has 
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led to the discovery of over 90 new halo GCs, extending to 
f^pmj ^ 140 kpc and 3D radii of >200 kpc (e.g, Huxor et al. 
2005; Huxor et al. 2008; Huxor et al 2013, di TulHo Zinn 
& Zinn 2013). A major step forward in understanding the 
formation of the outer halo GC system came from the real- 
isation that these objects preferentially lie on stellar streams 
and other debris features ( Mackev et al.ll20 1 Obi |20 1 3b . Monte 
Carlo simulations indicate a probability of < 1% that such 
alignments should happen by chance, leading to the conclu- 
sion that 80% of the M3 1 outer halo GCs have been accreted 
along w ith their ho st galaxies, confirming the idea put forward 
by Se arle & ZinnI (il978h for the Milky Way. In this Letter, we 
present the first results of a spectroscopic survey of these outer 
halo objects, focusing on their global kinematics. A detailed 
description of the data as well as a full analysis is deferred for 
a later publication (Veljanoski et al. in prep). 

2. THE DATA 

Spectra were acquired for 53 GCs spanning Rpmj ^20- 
130 kpc, of which 12(6) He beyond > 80(100) kpc. This 
sample is complete down to g 18.5, and 44 of the clus- 
ters had not previously been observed spectroscopically. The 
data were obtained over 15 nights during 2005-2010 using 
the ISIS spectrograph on the WHT 4.2m, and the RC spec- 
trograph on the KPNO 4m. ISIS has two detectors that in- 
dependently sample the blue and red spectral range. We 
used the R600B and R600R gratings to cover the wavelength 
range ^ 350-510 nm with a dispersion of 0.045 nm/pixel 
and ^ 750-920 nm with a dispersion of 0.079 nm/pixel re- 
spectively. With the RC spectrograph, we used the KPC007 
grating with a wavelength coverage of ^ 350-650 nm and a 
dispersion of 0.139 nm/pixel. Total integrations were 600- 
7200 seconds depending on the target brightness. The slit 
width was 1-2". The signal-to-noise per pixel was ss 7-30 
for most targets and 50-70 for the brightest objects. 

The data were reduced using standard IRAF0 procedures. 

IRAF is distributed by the National Optical Astronomy Observatories, 
which are operated by the Association of Universities for Research in Astron- 
omy, Inc., under cooperative agreement with the National Science Foundation 
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One-dimensional spectra were extracted with aperture radii 
of 2-2.5". Heliocentric radial velocities were derived using a 
chi-squared minimization technique between GC spectra and 
radial velocity template stars (Veljanoski et al. 2013, in prep). 
This is analogous to the standard cross-correlation method, 
and it produces similar results. The method has the advantage 
that it uses the uncertainties in both the target and template 
spectra, which helps to eliminate spurious features in the chi- 
squared function. 

The uncertainty in the radial velocity of each cluster is 
adopted to be the standard deviation of all the independent 
chi-squared minimizations between the cluster and multiple 
radial velocity standard stars. Furthermore, as we obtained 
two independent velocity measurements for the GCs observed 
with ISIS, these were combined to reduce the uncertainty in 
the final velocity value. The final median uncertainty of all 
53 GC velocity measurements is 12 km/s. Of the 9 GCs in 
our sample that had published velocities in the RBC, all agree 
to within one standard deviation. Four of these clusters were 
found to have more precise RBC velocities compared to our 
measurements and for these we adopt the RBC values in our 
subsequent analysis. 

As our GC sample spans a large extent on the sky, we con- 
verted our heliocentric radial velocities to the Galactocentric 
frame in order to remove any effects the Solar motion could 
have on the kinematics. This conversion was done usin g 
the relations presented in ICour teau & van den Bergh'(T999), 
with u pdat ed values for th e Solar motion from McMillan 
(1201 Ih and iSchonrich et all (12010). For the purpose of this 
study, w e take the M31 heliocentric ve locity to be -301 
+4 km/s (ICourteau & van den Berghll 19991) . which translates 
to a Galactocentric radial velocity of -109 +4 km/s. 

3. RESULTS 

3.1. The Global Velocity Map 

Figure [U shows the most recent metal-poor ([Fe/H] < -1 .4) 
red giant branch stellar density map from the Pan- Andromed a 
Archaeological Survey (PAndAS) ( McConnachie et al.l2009l) . 
Overlaid are the positions of the observed GCs, the colors of 
which correspond to their measured Galactocentric radial ve- 
locities. 

Three interesting groups of GCs are indicated on Figure 
[U The blue rectangle marks 4 GCs that project on the 
North-West stream. A contiguous velocity gradient is seen 
along this feature, with the most radially-distant cluster (Rproj 
«125 kpc) having a velocity of « -183 km/s and the in- 
nermost object (Rproj ~67 kpc) having « -380 km/s. Also 
marked (red contour) ar e the GCs t hat lie on stream D, a 
feature first identified by llbata et al.l (l2007h . There is no ap- 
parent velocity gradient in this case, but the velocities of the 
GCs observed in the northeastern part of the stream suggest 
two distinct kinematic groups. In particular, 4 GCs have 
a mean velocity and dispersion of -171 km/s and 50 km/s, 
while the remaining 3 GCs have mean velocity and disper- 
sion of 71 km/s and 18 km/s. Interestingly, the northeastern 
part of stream D has been previously identi fied to be a com- 
plicated region where it overlaps stream C ( Ibata et aP 120071 : 
[Richardson et al. 2011). Thus, it seems likely that the two 
kinematic GC groups we have identified are associated with 
these different streams, and are moving in the opposite sense 
around M31, as judged by their mean velocities. Radial ve- 
locities for stream stars in this overlapping region are not yet 
available, but the lone southern GC on stream D has a radial 



velocity that lies within ~ 2 km/s of stream stars tentatively 
identified in a nearby field (IChapman et al.ll2008h . Interest- 
ingly, the And-I dwarf, located on the south end of Stream D, 
has a Galactocentric velocity of ~'-190 km/s, similar to one 
of the kinematic GC groups. The And-IX dwarf, located on 
the north end of Stream D with velocity of ^-20km/s, appears 
uncorrelated with the GC groups. 

The yellow region on Figure [T]marks "association 2" which 
iMackev et al.l (|^10b) identified as a statistical overdensity of 
GCs not associated with any obvious underlying stellar de- 
bris feature. The spread of velocities indicates that not all of 
the GCs in this region can be members of a kinematically co- 
herent subgroup. However, measurements for the remaining 
^60% of the putative "association 2" GCs are required before 
the presence of a subgroup can be definitely ruled out. 

3.2. Rotation and Velocity Dispersion 

It has been known for some time that GCs in the inner 
regions of M31 rotate around the minor optical axis of the 
galaxy, in the same sense as the disk rotation. IPerrett et al.l 
(2002) measured a rotational amplitude of the GC system of 
- 140 km/s, while iLee etaP (2008) measured 190 km/s 
using a larger sample . Divi ding the GCs based on their met- 
alliticv. lDeason et aTl (1201 ll) found more pronounced rotation 
for the metal-rich ([Fe/H] >-l) subpopulation. Inspection of 
Figure[T]strongly suggests this rotation persists to larger radii, 
with GCs in the northeast having systematically higher veloc- 
ities than those in the southwest. Figure |2] shows the Galac- 
tocentric radial velocities, corrected for the systemic motion 
of M31, versus their projected distances along the M31 major 
axis. GCs belonging to the subgroups identified in Figure [T] 
are color-coded. The rotational signature appears to be a prop- 
erty of the bulk population of the outer halo GC sample and is 
not driven by one or two kinematically-coherent subgroups. 

To further inv estigate the rotational signature, we follow 
ICote et all ( 1200 l i) and fit the observed projected Galactocen- 
tric radial velocities Vp of the GCs with the function, 

vp(e) = v„.,,+Asin(e-eo) (i) 

where 6 is the projected position angle, measured east of 
north, of a cluster relative to M31 center, 6'() is the projected 
position angle of the GC system rotation axis, Vsys is the sys- 
temic velocity of the GC system and A is the amplitude of 
rotation. This approach assumes that the rotation axis of the 
GC system is perpendicular to the line of sight, and that the 
intrinsic angular velocity of the system is constant on spher- 
ical surfaces. The uncertainties are determin ed using the nu- 
merical bootstrapping technique (lEfronll982b and the derived 
rotational amplitudes are c orrected for the inclin ation of the 
M31 disk, taken to be 77.5°('Fergus on et al.ll2002h . 

We augment our radial velocities with those from the RBC 
and fit the GC sample as a whole as well as within and beyond 
30 kpc. This radius correspo nds to a clear break in the GC ra- 
dial number density profile (iHuxor et al.ll201 ll) and therefore 
provides a natural division between the "inner" and "outer" 
halo. For reference, the outer halo sample consists of the 48 
GCs presented here, to which we add a further 2 confirmed 
GCs from the RBC. 

The systemic velocity of the GC system was set to the M3 1 
Galactocentric systemic velocity. The results of this fitting are 
displayed in Table[T] The results remain unchanged when Vjvi 
in Equation |2] is left to vary as a free parameter, or when the 
mean velocity of the GC system is used. Given the position 
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Figure 1. The metal-poor stellar density map of M31 from PAndAS. Positions of the observed GCs are marked with colored dots which correspond to their 
Galactocentric radial velocities in units of km/s. Some GC subgroups are indicated (see text for details). The purple dashed circles con'espond to 30 and 130 kpc 
radii in projection. The Galactocentric velocity of M3 1 is - 109 +4 km/s. North is up and East is left. 



Table 1 

Derived Rotational Properties for M3 1 Halo GCs 





A 

[km/s] 


do 

[degrees] 


Velocity Dispersion 
[km/s] 


Ngc 


All GCs 


133+11 


124+4 


115 +5 


595 


Rpmj<30 kpc 


137+10 


124+4 


114+5 


545 


Rpro, >30 kpc 


79+19 


123 +27 


106+12 


50 



angle of the M31 major axis is 38°, the derived rotation axis 
is consistent with the minor axis of M31, and is essentially 
indistinguishable from the rotation axis of the inner halo GCs. 

Finally, we derived the rotation-corrected velocity disper- 
sion of the M31 GC sample using the biweight scale of 
|Beers et al. ( |199Ql) . Figure |3] shows the rotation-corrected 
Galactocentric radial velocities (corrected for the systemic 
motion of M31) versus Rproj of the M31 GCs. It reveals a de- 
creasing velocity dispersion with increasing distance from the 
M31 center, varying from ^122 km/s at 60 kpc to ^ 57 km/s 
at 120 kpc. For reference, the metal-poor field star velocity 



dispersion is ^98 km/s at 60 kpc (measured), a nd ~44 km/s 
at 120 kpc (extrapolated) (IChapman et al.ll2006l) . 

3.3. An M31 Mass Estimate 

Assuming the halo GC system is spherically symmetric, we 
can estimate the mass of M 3 1 by solving the Jeans equation 
jBinnev & Tremainell987 ). Because we have found evidence 
for rotation in the M3 1 halo GC system, the Jeans equation 
is separated into a rotating and a non-rotating component. 
The total mass is obtained by summing the mass supported 
by pressure Mp and the rotationally-supported mass M,-. The 
rotational component is determined via: 



R„ 



(2) 



where Rmax is the projected radius of the outermost GC in 
our sample, v„,ax is the rotational amplitude of the outer GC 
population and G is the gravitational constant. 

To determine Mp we use the solution of the non-rotating 
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system, it has the following form: 
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Figure 2. Galactocentric radial velocity, connected for the M3 1 systemic mo- 
tion, vs. projected radius along the M31 major axis. The black and colored 
squares mark the velocities of the GCs presented in this work, with different 
colors marking GC subgroups identified on Figure [T] RBC values are shown 
as small grey circles. The open squares correspond to the mean velocities 
of the GCs with Rp,oj>30 kpc calculated in 20 kpc bins. The .v error bars 
represent the bin size, while the y eiTor bars indicate the standard deviation 
of the mean. The green solid lines correspond to our measured amplitude for 
the outer GCs corrected for inclination. 
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Figure 3. Rotation-coiTected Galactocentric radial velocity, corrected for the 
M31 systemic motion, vs. Rpmj from the M31 center Symbols are as in Fig- 
ure[T] The open squares con'espond to the mean velocities of the halo clusters 
in 20 kpc bins. The x error bars mark the bin size, while the y error bars rep- 
resent the rotation-coiTected velocity dispersion. The dotted line marks the 
30 kpc radius. 



Jeans equation proposed by IE vans et alj (12003b . known as the 
tracer mass estimator (TME): 



M„ 



GN 



(3) 



where Ri is the projected radius from the center of M3 1 for 
a given cluster, v, is the radial velocity of the GC with the 
rotational component removed and is the total number of 
clusters in our sample. The constant C depends on the shape 
of the potential, the radial distribution of the tracer objects 
and the anisotropy in the system. For a spherical, isotropic 



C = 



4(a-l-7) 4-a-7 l-{A,Jr„ud 
TT 3-7 l-(nn/rout)^ ■ 



(4) 



In equation (01, rjn and rout correspond to the smallest and 
largest 3D radii of the halo GCs respectively. In this case, 
we assume r-m = 30 kpc, while rout is set to 200 kpc assum - 
ing MGCl is the most remote GC (iMackev et alj I2010al) . 
The constant a is related to the underlying gravitational field, 
which is assumed to be scale-free, at least between rin and rout- 
For an isothermal halo potential where the system has a flat ro- 
tation curve at large radii, a is zero. If we assume a NF W dark 
matter profile (iNavarro et a l. 1996 ), a « 0.55 ( Watkin s et alJ 
I2OIQ1) . The 7 parameter is the slope of the GC volume den- 
sity distribution. We calculate this using the surface density 
distribution of all 84 GCs that have projected radii larger than 
30 kpc (Huxor et al. 2013, in prep; Mackey et al. 2013, in 
prep) and find 7 ^ 3.34. It is also important to note that even 
though the TME uses a GC sample in a shell around the cen- 
ter of M3 1 , it calculates the total mass enclosed by the furthest 
cluster in that sample. 

Using the above method, and setting a to zero, we calcu- 
late the total mass of M31 tobeMM3i = 1.5 ±0.2 x IO^^Mq, 
where the pressure component is Mp = 1.3 ±0.2 x IO'^Mq 
and the rotation contribution is = 2 ± 1 x lO^M©. As- 
suming a = 0.55, we find the M31 mass to be Mm3i = 1.2 ± 
0.2 X lO'^M©, where Mp = 1.0 ±0.2 x lO'^M© and M,- = 
2 ± 1 X 10"Mq. For reference, applying the TME in a sin- 
gle step (ignoring rotation), gives Mm3i = 1 .8 ± 0.2 x IO'^Mq 
andMM3i = 1.3 ±0.2 x IO'^Mq for a values of and 0.55 re- 
spectively. The quoted errors incorporate the statistical uncer- 
tainties only, and in reality they are much larger In our mass 
calculations, we assume isotropic orbits, a steady state for our 
tracer population and a power-law form for the potential. This 
is the simplest approach we can take, although the presence 
of the substructure in the spatial distribution of the outer halo 
GCs suggests that it may not be correct. Nonetheless, stud- 
ies suggest the presence of substructure in the t racer popula- 
tion will bias results only at the 20% level (e.g.. lYencho et al] 
1^06; Deason et al. 2012). To explicitly test the assumption 
of steady state, we recalculate the M3 1 mass excluding GCs 
which lie along the North-West stream and stream D. We find 
the total mass of M31 decreases by 0.3 x lO'^M© for both 
values of a, with the formal statistical errors remaining un- 
changed. 

4. DISCUSSION 

Analysis of the radial velocities of M31 outer halo GCs 
strongly supports our earlier finding that many of these ob- 
jects have been accreted (Mackey et al. 2010b). Clear kine- 
matic correlations are seen amongst subgroups of GCs which 
lie on top of stellar debris features, and, in the case of the 
North-West stream, an unambiguous signature of radial in- 
fall is observed. Interestingly, recent work has also began 
to detect phase-space substruc ture in large samples of GCs 
around distant galaxies (e.g. S trader et al.ll201 It IBlom et aP 
2012; Romanowsky et al. 2OI2F" However it is only within 
the Local Group that we can attempt the obvious next step 
of comparing GC velocities with those of underlying stream 
stars. 

A surprising discovery is the high degree of rotation in 
the M31 outer halo GC population, which is in the same 
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sense as the inner halo population as well as the main stel- 
lar disk. While it is common to find GCs rapidly rotating in 
inner regions of galaxies, strong rotation beyond a few tens 
of kpc seems to be a rare occurrence, at least amongst early 
type galaxies (e .g. Woodlev et al.' ' 20101: ISfa-ader et all [201 It 
iBlom et all 12012 : Potaetal. 2013). It is nahiral to specu- 
late on how such coherent motion could arise if the outer 
GC population were largely accreted from numerous dwarf 
galaxy hosts. One possibility is that the donor dwarf galax- 
ies have been accreted from a few preferred directions on 
the sky and hence have aligned angular momentajjis seen in 
some recent cosmo logical simulations (iLibeskind et^ni2005L 
l20TltlLovelletani2011.) . Such a scenario has also been sug- 
gested as the origin of the planar alignments of dwarf galaxies 
seen in bot h the Milky Way and M31 (e.gj Metz et alj|2007l: 
llbata et a"l] r2013). Curiously, the plane of dwarf galaxies re- 
ported in M3 1 rotates in the same sense as the outer halo GC 
population, although the rotation axis of that plane appears to 
be inclined by ^ 45deg to the minor axis. It is also possi- 
ble that the bulk of the M31 outer halo GC population was 
accreted in a single event involving a moderate mass satel- 
lite. Support for this idea could come from the M31 thick 
disk which rotates in the sa me sense, albeit som ewhat faster, 
than the halo GC population dColUns et al.ll201 lb . However, a 
rather massive satellite would be required to bring in a popu- 
lation of several tens of GCs, raising the question of how the 
M31 disk could survive such an encounter Furthermore, it 
would seem difficult to explain the spatial correlation between 
outer halo GCs and the numerous tidal streams in this case. In 
a different scenario, numerical modelling (Be kki 201 0) sug- 
gests that a past major merger between M3 1 and another disk 
galaxy could give rise to the rapid rotation of the resulting GC 
system, including the rotation in the halo population. 

Despite being our closest massive galaxy, it is remarkable 
that we are still unable to measure the total mass of M3 1 with 
good precision. Indeed, there is still deba te as to whethe r 
M3 1 or the Milky W ay is more massive (Watkins et aLllIoTO ). 
IE vans et aTl ( l2003h used GC kinematics to find a M3 1 mass 
of 1.2 X l O'^Mfn out to a deprojected 3D distance of ^^100 
kpc, while iGalleti et al. (2006) and Lee et al. (2008) used ex- 
panded samples to calculate masses of 1.9-2.4 x IO'^Mq 
within the same radial range. Several authors have attempted 
to determine the M3 1 mass via the motions of its dwar f satel- 
lite galaxies (e.g. ICdte et al.ll2000l: Evans et al."2003' ). The 
most recent such measurement is presented by Watkin s~et al.l 
(I201QI) who estimated 1.4 ±0.4 x IO'^Mq using 23 satel- 
lites out to 300 kpc, assuming isotropy. Our estimate of 
Mm3i = 1.2- 1.5 ±0.2 X lO'^M© within a 3D radius of 200 
kpc agrees well with this value but suffers from the simi- 
lar syst ematic uncertainties due to model assumptions. No- 
tably, van der Marel et akl (120121) use the velocity vector of 
M3 1 with respect to the Milky Way in combination with the 
timing argument to derive a total mass for the Local Group 
Mlc = (4.39 ± 1.63) X IO'^Mq, which would push the M31 
mass higher than any estimate thus far using dynamical trac- 
ers. 
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